AD-A122  255  AN  APPROXIMATION  FOR  G(S)F(S))  VIA  THE  SECOND  MEAN 
VALUE  THEOREM  OF  THE  I..(U)  ARMY  MISSILE  COMMAND 
REDSTONE  ARSENAL  AL  SYSTEMS  SIMULATION  A..  RE  DICKSON 
UNCLASSIFIED  NOV  82  ORSMI/RD-82- 16-TR  S8I-A0-E95O  322  F/G  12/1 


*-  :  <  t 


Q 

<XT 


?  •  v 


TECHNICAL  REPORT  RD-82-16 


AN  APPROXIMATION  FOR  ZtG(s)F(s)J  VIA  THE  SECOND 
MEAN  VALUE  THEOREM  OF  THE  INTEGRAL  CALCULUS: 

THE  DIGITAL  SIMULATION  OF  CONTINUOUS  LINEAR 
SYSTEMS  VIA  Z-TRANSFORMS 


Richard  E.  Dickson 

Systems  Simulation  and  Development  Directorate 
US  Army  Missile  Laboratory 


November  1982 


dDJ'AIM  MDSSOILE  COMMAND 

Rec l&tone  XKr&&r*m.l ,  Afabama  35898 


>- 

CL. 

C 

<L 


L,t 


K  * 

r*. 


Approved  for  public  release;  distribution  unlimited. 

DTtC 

^1-ET.CTE 


lV'"0EC  10  1982 


82  12  08  014 


aw  FORM  1081,  1  JUL  79  PREVIOUS  EDITION  IS  OBSOLETE 


Unc  1  ass i  t  it'd 


SECURITY  CLASSIFICATION  OF  THIS  PACE  'WNl.n  Pete  EnM.mil' 


REPORT  DOCUMENTATION  PAGE 


1  REPORT  NUMBER  2  GOVT  ACCESSION  I 

RD-82-16 

«  Title  f«nd  Subtitle) 

AN  APPROXIMATION  FOR  Z[G(s)F(s)]  VIA  THE  SECOND 
MEAN  VALUE  THEOREM  OF  THE  INTEGRAL  CALCULUS: 

THE  DIGITAL  SIMULATION  OF  CONTINUOUS  LINEAR 


Sire  .  READ  INSTRUCTIONS 

MU1-  _ _ before  completing  form 

2  GOVT  ACCESSION  NO.!  3  RECIPIENT'S  CATALOG  NUMBER 


5.  TYPE  OF  REPORT  &  PERIOD  COVERED 


Technical  Report 


6.  PERFORMING  ORG.  REPORT  NUMBER 


r7.  AUTHOR(A) 

Richard  E.  Dickson 


8.  CONTRACT  OR  GRANT  NUMBERCM 


9-  PERFORMING  ORGANIZATION  NAME  AND  ADDRESS 

Commander,  US  Army  Missile  Command 
ATTN :  DRSMI-RD 


10.  PROGRAM  ELEMENT.  PROJECT.  TASK 
AREA  &  WORK  UNIT  NUMBERS 


11.  CONTROLLING  OFFICE  NAME  ANO  ADDRESS 

Commander,  US  Army  Missile  Command 
ATTN:  DRSMI-RPT 


12.  REPORT  DATE 

November  1982 

13.  number  of  pages 


U.  MONITORING  AGENCY  NAME  6  ADDRESS  (It  dltlerent  I  ram  Controlling  O I  lice)  IS.  SECURITY  CLASS,  (of  title  report) 

Unclassif ied 

I5».  DECLASSIFICATION/  DOWNGRADING- 
SCHEDULE 

US.  DISTRIBUTION  STATEMENT  (ol  thle  Report)  - - __ - _____ — 

Accession  For 

Approved  for  public  release;  distribution  unlimited.  NTIS  GRAM  ~~ 

DT1C  TAB  q* 

Unannounced  Q 

Justification _ _ 


17.  DISTRIBUTION  STATEMENT  (ol  the  mbetrect  entered  In  Block  20,  II  dltlerent  tram  RepoA 


18.  SUPPLEMENTARY  NOTES 


1  °/-v 


By - _ - 

Diet rl but 1 on/ 
Availability  Codes 


)  Dist  Special 


19.  KEY  WORDS  (Continue  on  rtr«r«#  mi  dm  it  nmc+mmary  and  identify  by  btock  number) 

Z-transforms 
Digital  simulation 
Linear  Systems 

20.  ABSTRACT  (C oottmim  «n  rmrmrmm  midm  ft  nmeammary  and  idmtUlfy  by  block  numbar) 

This  report  presents  the  approximation  for  Z[G(s)F(s)],  which  is  based 
upon  the  second  mean  value  theorem  of  the  integral  calculus. 


DD  .KT*  1473  corn  on  of  t  nov  ss  is  obsolete 


Unclassified 

SECURITY  CLASSIFICATION  OF  THIS  PAGE  fiwi  Dale  Entered) 


AN  APPROXIMATION  FOR  Z[G(s)F(s)]  VIA  THE  SECOND  MEAN  VALUE  THEOREM  OF  THE  INTEGRAL  CALCULUS: 
THE  DIGITAL  SIMULATION  OF  CONTINUOUS  LINEAR  SYSTEMS  VIA  Z-TRANSFORMS 


Richard  E.  Dickson 

Systems  Simulation  and  Development  Directorate 
US  Army  Missile  Laboratory 
US  Army  Missile  Command 
Redstone  Arsenal,  AL  35898 


Summary 

The  problem  with  applying  z-transforms  to  the  digi¬ 
tal  simulation  of  continuous  linear  systems  is  that  the 
Raggazzini-Zadeh  identity  only  applies  to  the  initial 
conditions.  This  difficulty  may  be  overcome  by  approx¬ 
imating  Z[G(s)F(s)),  for  example,  Halijak's  "Trapezoidal 
Convolution"1 and  its  extensions.3  These  approxima¬ 
tions  did  not  yield  perfect  response  for  a  unit  step 
into  a  single  pole  filter  unless  adjusted  by  a  residue." 

The  approximation  for  Z[G{s)F(s)]  presented  here, 
which  is  based  upon  the  second  mean  value  theorem  of 
the  integral  calculus,  has  the  desired  response  for  a 
unit  step  into  a  single  pole  filter  and  does  not  require 
adjustment. 


Some  Prel iminaries 


The  "sifting"  property  of  the  Dirac  delta,  equation  (4) 
yields 

Z[G(z)F(s)l  =  t  I  g(kT)z‘kf(nT  -  kT)z'(n'k)  (9) 
n=0  k=0 

where  the  summation  over  k  is  discrete  convolution. 
From  the  Cauchy  product  of  power  series 

z  G(z)F(s)]=^  £q  g(kT)z"kj  ^  £  f(nT)z‘nj  (10) 

and  from  equation  (7) 

Z[G(z)F(s)]  =  G(z)F(z)  (11) 

the  Raggazzini-Zadeh  identity. 


The  Laplace  transform  is  defined  to  be 


F(s)  =  /  f(t)  u(t)  e'st  dt 

— oo 


(1) 


where  u(t)  is  the  Heaviside  unit  step. 


u(t) 


I  0,t  <  0 
1  l,t  >  0 


(2) 


The  z-transform  will  be  defined  as  a  "discrete"  Laplace 
transform 


F(z)  =  /  1  f(t)  5(nT  -  t)  u(t)  e*st  dt  (3) 

n=-°° 

and  from  the  "sifting"  property  of  the  Dirac  delta  dis¬ 
tribution, 


Z[G(s)F(s)l  Approximation 

It  is  well  known  that,  in  general,  Z[G(s)F(s}] 
i  T  G(z)F(z)  and  therein  lies  the  difficulty  in  apply¬ 
ing  z-transforms  to  the  digital  simulation  of  continu¬ 
ous  systems.  A  very  useful  relationship  would  be  the 
solution  of 


z[s(s)F(s)]  «  I  z'n  u(nT) 

n- -  oo 


X  f  g ( t)  u ( t)  f(nT  -  t)  u(nT  -  t)  dt 


If  both  G(s)  and  F(s)  are  known  a  priori ,  there  is 
no  difficulty;  but,  in  digital  simulation,  though  the 
plant,  G(s),  is  known  a  priori,  the  input,  f(nT),  is 
not. 


f(nT)  =  /  f(t)  5(nT  -  t)  dt 

-CO 

equation  (3)  becomes 

F(z)  =  I  f(nT)  u(nT)  e'snT 
n--°° 

and  letting 
one  has 


z  =  e 


(4) 

(5) 

(6) 
(7) 


By  the  second  mean  value  theorem  of  the  integral 
calculus.  Bonnet's  theorem,  one  may  write 

”  n  1  i  kT+n^T 

Z[G(s)F(s)]  =  l  t  ,f(nT  -  kT)  /  g(t)  dt 

n=0  k=0  I  kT 


+  f(nT  -  kT  -  T) 


kT+T 


/  9 

kT+nJ 


g(t)  dt  z 
kT  \ 


(13) 

-n 


F(z)  =  1  f(nT)  z‘n 

n=0 

the  usual  definition2  of  the  z-transform. 

The  Raggazzini-Zadeh  Identity  may  be  readily  deduced 
from  the  convolution  properties  of  the  Laplace  transform 


In  general  n^  would  be  different  for  each  k 
Assume  that  n  is  the  same  for  each  k  ,  that  is. 


nk  -  n 


(14) 


Equation  (13)  is  almost,  but  not  quite,  the  form 
of  the  Cauchy  product  of  power  series.  To  obtain  the 


Z(G(z)  F(s)]  *  £ 

n*-« 


u(nT)  f  £  g(t)  6(kT  -  t)  u(t)  f(nT  -  t)  u(nT  -t)  dt 

CT>  |(S— CD 


(8) 


1 


proper  form 


“  (  n-1  kJ+nT  n  f  -n 

Z[G(s)F(s)]  =  I  ,  I  f (nT  -  kT)  f  g(t)dt  +  I  f(nT  -  kT)  /  g(t)  dt  -  z 

n=0  |  k=0  kT  k=l  kT+(n-l)T  \ 


Z  [G  ( s )  F  ( s )]  =  1  \  t  f (nT  -  kT)z*(n'k)  /'g(t)dt  z‘k  -  f(0)/"g(t) 

n=0  | k=0  kT  nF 


kT +nT 


nT+nT 


dt  z" 


+  f  f (nT  -  kT)z'(n'k)  /  Tg(t)dt  z'k  -  f(nT)  z'n  /  g(t)dt| 
k=0  k T + ( n - 1 ) T  (n-l)T  \ 


and  finally 


nT+nT 


Z  [G(s)F(s)]  ®  [F(z)  -  f ( o )]  £  f  g(t)  dt  z  +  F(z) 

n=0  nT 


nT 


I  /  g(t)dt  z'n  -  /g(t)dt 
n=0  nT+(n-l)T  (n-l)T  J 

It  is  interesting  to  compare  this  result  with  Halijak's  "Trapezoidal  Convolution:"1’2 

Z[G(s)F(s)]  =  |-  { [F(z)  -  f(o)]  G(z )  +  F(z)  [G(z)  -  g(o)]J 


which  may  be  derived  from  the  first  mean  value  theorem 
of  the  integral  calculus.’ 


Two  Plants  of  Interest 

Consider  the  case  where  the  plant  is  a  single  inte¬ 
grator;  that  is: 

G(s)  =  1/s  (19a) 

g(t)  ■  1  (19b) 

Substituting  into  equation  (18)  one  has 

«°  nT+nT 

Z[F(s)/s]  =  [F  ( z )  -  f  { o )]  If  dt  z'n 

n=0  nT 


that  is 


and 


f(t)  *  1 


F(z)  =  1/(1  -  z'1) 

equation  (23)  becomes 

d-z'T  i  -  z*‘ 

2  [(l/s)/s] 


+  F(z) 

and  integrating 


«  nT  „  0 

£  /  dt  z  -  f  dt 
n=0  nT+(n -1  )T  (n-l)T  J 


(20) 


Z[ F(s)/s]  =  0 (z)  -  f(o)]  t  z-n 

n-0 


(1  -  z’;)2 

which  is  exact  for  any  r,  . 

For  F(s)  =  1/s- 

that  is  f(t)  =  t 

and  F(z)  =  Tz/(  1  -  z'1  )2 

equation  (23)  becomes 


+  F(z) 


I  (l-n)T  z-"  -  (l-n)T 
n=0 


Z[(l/s2)/s] 


..  [r!  t  U  ~  "Oz'jlT  /  T  z 


(21) 


1  -  z 


(CT) 


which  may  be  written  after  summing 
Z[F(s)/s]  =  [F ( z )  -  f(o)]  [nT/(l-z"! )] 

+  F(z)  [(l-n)T/(l-z“)  -  (l-n)T] 
and  combining  like  terms 

Z[F(S)/S]  -  [n  ±  (l-n)z-]  T  F(zl  -,nT.f(0) 

1  -  z 

For  a  unit  step  input, 

F(s)  =  1/s 


(22) 


(23) 


(24a) 


z[(i/s?)/s]  = 

(I  -  z'1)5 

which  for  n  *  1/2  is 

z[(i/s2)/s]  =  JIiZLLljI!! 

2(1  -  z->)3 

which  is  exact. 

A  more  interesting  plant  is  the  single  pole 
G(s)  =  a/(s  +  a) 

that  Is 


(15) 

(16) 

(17) 

(18) 

(24b) 

(24c) 

(25) 

(26) 

(27a) 

(27b) 

(27c) 

(28) 

(29) 

1 30 ) 

filter: 

(31a) 


9(t)  =  a  e'aT 

In  this  case  equation  (17)  becomes 


(31b)  filter  when  the  time  constant  of  the  input  and  the  fil¬ 
ter  are  the  same1' 


_  _  ~  1 1  1  '  '  I  I  , 

Z[a  F(s)/(s  +  a)]  =  [F(z)  -  f(o)]  I  f  ae‘at  dt  z'n  +  F(z) 

n=0  nT 

and  after  performing  the  indicated  integration  and  sums 


Z[a  F(s)/(s  +  a)]  =  [F(z)  -  f (o)] /-1-  +  F(z) 

(l  -  e'at  z'1) 


“  nT 

I  /  ae’atdt  z’n  - 
_n=0  nT+(--l)T 

Cf,J(e-a("-1)T-  Ik81 


atdt  z'n  -  f  ae*at  dt  (32) 
1 )T  (n-l)T  J 


1  -  e"at  z- 


In  this  case,  "Trapezoidal  Convolution,"  equation  (18)  would  yield 

Z[a  F(s)/(s  +  a)]  =  (|^)J[f(z)  -  f(o)]  /  (l-e"dl  z*1)  +  F v. z )  [e'aT  z‘!/(l  -  e"aT  z'1)]! 

and  with  residue  adjustment 

Z [a  F(s)/(s  +  a)]  =  tanh(^)  { [ F ( z )  -  f(o)]/(l  -  e~aT  z*1)  +  F(z)  [e~aT  z‘‘/(l-e'aT  z")j 

For  a  unit  step  input,  F(s)  =  1/s  (24a) 

7r  1  /  a  \  aT  e"a  z’1 

that  is  f(t)  =  1  (24b)  +  a  Is  +  a/J  '  (j  .  e-aT  z-TjT 

and  F(z)  =  l/(l  -  z"1)  (24c)  In  this  case,  equation  (17)  yields 

equation  (33)  becomes  r  i  /  a  Cl  e'a!,T(l  -  e'aT)z‘‘ 

i  n  i+) l  -  r^,  •  .i  \  ’  o  •  .-T 


[-1  (r*h)]  ■  [rh-  ‘  >]  (^r,.) 


z[rh(rb)  ■ 

this  case,  equati 

‘[rbfciriM 


aT  e‘aT  z“ 

(r."e-aTz 


F(z)  =  1/(1  -  z  ')  (24c)  In  this  case,  equation  (17)  yields 

a  \l  e’ar'T (l  -  e‘aT)z’1 


“  1  (1  -  e  Vz 

(1  -  e’aT  z-r 


+  (l 

1  N  {"(e"anT  -  e-aT)z-‘l 

-  ,  .  ,-T  J 

(36) 

collecting  terms 

z[s  (.  .*.)]■ 

(l  - 

l  *37  \ 

(l  -  z-')(l  -  e'aT  z'') 

W'  J 

and  setting  this  numerator  equal  to  the  exact  numerator 
and  solving  for  n  ,  one  has 


.  =1  IT  e 


n  ~  aT  tn 


For  aT  <<  1 

n  ~  1/2 


which  is  exact  for  any  o  ! 

Without  residue  adjustment  "trapezoidal 
convolution"  yields1* 


l  _aj\  -l  cer  W1ln  non-zero  ini 

zTi  /  a  VI  -  /aT\  (1  +  e  )z _ differential  equation 

Ls  'S  +  a/J  '  \2  /(  .  7%  .  e-aT7T 

V  Z  /'  Z  ’  x(t)  +  a  x(t)  =  a  f  (t) 

With  residue  adjustment  one  has 

n  /  a  M  /aT\  (l  +  e'aT)  z"1  ,  and  transform  to  the  Laplace  domain  using 

%(rrr)J  ■  •»»  (r)(l .  ^  ( e-ii (39i> 

zri/_i_jl  .  ii  -  .-I.-  „9bl  3  ■ 

(i  -  rVi  - r1) 

'  ''  s  X ( s )  -  x (0)  + 

In  the  sample  data  domain  a  second  order  Runge-Kutta 

integrator  would  yield:  solving  for  X(s)  and 

zfl  (  a  S\_  _ Q  -  0  -  »T  ♦  a2T Vgfl_? _ _  (40) 

Is  Is  +  aJT  (1  .  z  >)[i  -  (1  -  aT  +  a2T2/2)  z‘‘]  x(z»  ,  z  [aF(. 


Recurrence  for  A  Single  Pole  Filter 

To  develop  a  recurrence  for  a  single  real  pole  fil¬ 
ter  with  non-zero  initial  conditions,  start  with  the 


[<<”>(,)] 


snF(S)  -  "iV-'-’fio! 
s=o 


s  X(s)  -  x ( 0 )  +  a  X ( s )  =  a  F(s) 
solving  for  X(s)  and  then  taking  the  z-transform. 


X(z)  =  Z 


Note  the  (2/0)  Pade"  approximation  for  the  exponential. 

A  unit  step  input  and  a  Tustin's  substitution2  for 
the  plant  would  lead  to 


tm] 


Since  the  initial  condition,  x(o)  ,  is  a  constant 
in  the  Laplace  domain  (an  impulse  in  the  time  domain), 
the  Raggazzini-Zadeh  identity.  Equation  (11),  would 


;ri/  a  ft  -  (1  -  aT/2 )/ ( 1  ±  aT/2)J  z~  -  (41)  apply.  The  forcing  function,  F(s)  ,  is  another  matter, 

LsVs  +  a/J  (j  .  2-l)[i  -  (1  -  aT/2)z'J(l  +  aT/2)J  and  requires  approximation.  Using  equation  (33) 

Note  the  (1/1)  Pade"  approximations  for  the  exponential.  ^  _z-^-aT)x  (z )  -  [(l-e"anT)  +  (e'anT  -  e‘aT)z'*]  F(z) 

For  the  large  time  step,  T,  the  approximation  pre-  *- 

sented  here  will  out  perform  "Trapezoidal  Convolution,1'  _  (ie-anT\  ,  ,  >  (50) 

Tustin's  substitution  method  and  a  second  order  Runge-  '  '  ' 

Kutta  Integrator.  an(j  substituting  equation  (7)  and  equating  coefficients 

"Trapezoidal  Convolution”1’2  yields  exact  results  of  like  powers  of  "z" 
for  a  damped  exponential  input  into  a  single  real  pole  x(0)  =  x(0)  (51a) 


Conclusions 


x(nT)  *  e"aT  x(nT-T)  +  (l  -  e'anT)  f(nT) 

+  (e"anT  -  e'aT)f(nT  -  T),  n  >  0  (51b) 

For  a  second  order  damped  filter  with  non-zero 
initial  conditions4 

F(s)  +  (s  +  2o)  x ( 0 )  +  x(0) 

X  ( s )  =  -9 - (52) 

S2  +  2os  +  i ii2 
o 

where 

0  <  <  i 

U) 

0 

a  partial  fraction  expansion  yields 

*  2w  [s  +  a  +  iu  +  s  +  a  -  iw  ]  * 

[u2  F(s)  +  x(0)  +  x (0)] 

+  KrrTTTa, +  ttWiz]  *(0)  <53> 

where 

<*>  =  (<J2  -  02)1/2  (54) 

Let 

X ( s )  x  ill!  (55) 

where 

i  fu2  F(s)  +  a  x (0 )  +  x(0)l/ij  +  x(0) 

U<s>  =  — - - - -  <56> 

and 

-ifu2  F(s)  +  o  x(0)  +  x(0)l/u  +  x(0) 

V(s)  =  - -3  +-q--iS-  -J -  (57) 

Since  V ( s)  is  the  complex  conjugate  of  U(s)  , 
only  u(nT)  need  be  found.  The  problem  becomes  that 
of  finding  the  recurrence  for  a  single  complex  pole 
filter.  In  this  case,  equation  (51  a,b)  would  become: 

u(0)  -  x(0)  +  i[ax(0)  +  x(0)]/w  (58a) 

u(nT) 

where 

e-(o+iu)T  =  g-aT  ^CQS  _  .  s^n  (59) 

and 

x(nT)  =  Real  u(nT)  (60a) 

x(nT)  =  to  Imaginary  [u(nT)]-  a  Real[u(nT)]  (60b) 


The  approximation  proposed  here  is  exact  for  a 
unit  step  input  into  a  single  pole  filter.  The  pole 
may  be  either  real,  imaginary,  or  complex.  The 
response  to  a  unit  step  is  unaffected  by  the  choice  of 
n  ;  and  n  may  be  chosen  to  adjust  the  response  for 
some  other  input  of  interest,  though  a  value  of  one 
half  is  generally  appropriate. 

Higher  order  problems  should  first  be  “sirnpl ified" 
by  a  partial  fraction  expansion,  though  the  initial 
conditions  should  be  incorporated  before  the  expansion. 
The  approach  lends  itself  to  parallel  processing. 
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=  e-(a+iu}T  u(nT  -  T)  +  iw(J)2  [l-e'n(a+iw)T]  f(nT)  ♦  iu(jjf' 


n(r+iw)T  g-(o+iu))T 


]]f(nT  -  T) 


(58b) 
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